#include <stdio.h>
#include <stdlib.h>
#include <stdint.h>
#include <sys/time.h>

#define TRUE 1
#define FALSE 0

#define DEBUG FALSE

long calcRuntime(struct timeval *begin, struct timeval *end);

double drand();

void initMatrix(double* matr,int dim,int random);

void mult(double *a,double *b,double *c, int n);

void plot(double* mat,int dim);

int main(int argc,char** argv){
	
	double *a;
	double *b;
	double *c;

	struct timeval start, end;
	long runtime = 0.0;
	float runtime_sec = 0.0;
	float flops = 0;
	
	unsigned int size = 0;

	if(argc == 2){
		size = atoi(argv[1]);
	}else{
		size = 2048;
	}

	#if DEBUG == TRUE
		printf("Start calculation with size == %d\n", size);
	#endif

	a = (double*)malloc(size*size*sizeof(double));
	b = (double*)malloc(size*size*sizeof(double));
	c = (double*)malloc(size*size*sizeof(double));
	
	initMatrix(a,size,TRUE);
	initMatrix(b,size,TRUE);
	initMatrix(c,size,FALSE);

	#if DEBUG == TRUE
	printf("Matrix a:\n");
		plot(a,size);
	printf("Matrix b:\n");
		plot(b,size);
	#endif
	
	gettimeofday(&start,NULL);

	mult(a,b,c,size);

	gettimeofday(&end,NULL);
		
	#if DEBUG == TRUE
	printf("Matrix c:\n");
		plot(c,size);
	#endif

	runtime = calcRuntime(&start,&end);
	runtime_sec = (float)runtime/1e6;
	flops = (2*size*size)/(runtime_sec);
	
	
	printf("#runtime[us] runtime[sec] flops size\n");
	printf("%ld %f %f %d\n",runtime,runtime_sec,flops,size);

	free(a);
	free(b);
	free(c);

	return 0;
}/*END OF MAIN*/

/*Function: calcRuntime*/
long calcRuntime(struct timeval *begin, struct timeval *end){
	struct timeval diff;
	diff.tv_sec = end->tv_sec - begin->tv_sec;

	if(begin->tv_usec > end->tv_usec){
		diff.tv_usec = (1e6- end->tv_usec) + begin->tv_usec;
	}else{
		diff.tv_usec = end->tv_usec - begin->tv_usec;
	}

	return (diff.tv_sec * 1e6 + diff.tv_usec);
}

/*Function: drand()*/
double drand(){
	return rand() / (RAND_MAX + 1.);
}

/*Function: initMatrix()*/
void initMatrix(double* matr,int dim,int random){
	int i;
	int j;
	
	if(random == TRUE){
		for(i = 0; i <dim; i++){
			for(j = 0; j < dim; j++){
				*(matr+j+i*dim)=drand();
			}
		}
	}else{
		for(i = 0; i <dim; i++){
			for(j = 0; j < dim; j++){
				*(matr+j+i*dim)=0.0;
			}
		}
	}
}

/*Function: mult()*/
void mult(double *a,double *b,double *c, int n){
	int i; //Colums
	int j; //Rows
	int k; 

	//Berechnung 
	for(i=0;i<n;i++){
		for(j=0;j<n;j++){
			for(k=0;k<n;k++){
				*(c+j+i*n)+=(*(a+k+i*n))*(*(b+j+k*n));
			}
		}
	}
}

/*Function: plot()*/
void plot(double* mat,int dim){
	int i;
	int j;
	
	for(i=0;i<dim;i++){
		for(j=0;j<dim;j++){
			printf("%lf ",*(mat+j+i*dim));
		}
		printf("\n");
	}	
	
}
